Inverse problems: Can we obtain more? Quantum dynamics and the "^He case 



O 
(N 



E. Vitali, M. Rossi, L. Reatto, and D.E. Galli 

Dipartimento di Fisica, Universita degli Studi di Milano, via Celoria 16, 20133 Milano, 

(Dated: March 17, 2010) 



Italy 



We introduce a new strategy, called Genetic Inversion via Falsification of Theories (GIFT), to 
face inverse problems and we apply it to the extraction of information about real time dynamics of 
a many-body quantum system from noisy imaginary time correlation functions /(t) computed via 
Quantum Monte Carlo (QMC). Production and falsification of model spectral functions s{ijj) are 
obtained via a survival-to-compatibility with /(r) evolutionary process, based on Genetic Algo- 
rithms. Statistical uncertainty in /(r) is promoted to be an asset via a sampling of equivalent /(r) 
within the noise, which give rise to independent evolutionary processes. We have studied bulk *He 
at T = K; for the first time we recover from exact QMC simulations sharp quasi-particle excita- 
tions with spectral functions displaying also the multiphonon branch. As further applications of the 
method we have studied the impuriton branch of one ^He atom in liquid ''He and the vacancy- wave 
excitations in hep solid *He finding a novel roton like feature. 
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I. INTRODUCTION 

Since the earliest days of research in Physics, inverse 
problems have always provided challenges in a huge vari- 
ety of physical or even more generally scientific studies^^S.. 
At the most general level, an inverse problem consists in 
tising experimental observations in order to infer the val- 
ues of the parameters of a theory describing a natural 
phenomenon. At a first glance, one immediately con- 
vinces himself that such an inverse procedure in most 
realistic situations is unavoidably ill-posed, since any set 
of observations is limited and noisy, thus ruling out the 
possibility of finding out one and only one theory whose 
predictions fit such data. In more abstract terms, the 
problem of inverting a map connecting a space of theo- 
ries, say S, and a space of data, say T), is mathematically 
ill-posed; the reason for this lies in the inescapable limi- 
tations related to any realistic way of performing obser- 
vations, making impossible to have access to a sufficiently 
exhaustive and accurate knowledge of the space V. 

Some important questions thus naturally arise: what 
can we really learn when facing an inverse problem? 
What do we mean when we speak about a solution'! In 
our opinion, a key proposition is brilliantly put forward 
in ReL^i, following Popper**: observations may be used 
only to falsify a theory. In other words, we cannot ex- 
pect to find out a recipe that will allow to deduce all 
what is needed to build up a unique theory from a given 
set of observed data. Nevertheless, provided that we are 
able to construct a suitable parametrization of the ab- 
stract space of theories S, we may use observations to 
provide a "falsification test", aimed to exclude the the- 
ories whose predictions fail to fit the data. In this way 
we will be able to collect a (maybe very big) class of 
theories, which have been not falsified by the measured 
data. In our opinion, the best way to achieve this is to 
fully exploit all the information related to the observa- 
tions: that is, since any set of experimental data appears 
together with statistical uncertainties evaluated starting 



from suitable measurements, any set of data compatible 
with the original one has to take part to the falsification 
test, in order to suppress the possibility of unphysical 
effects arising from statistical fluctuations. 

Once we remain with a set of equivalent "solutions" 
"survived" to the falsification test, depending on the 
mathematical details of the space S, a natural idea ap- 
pears to be that of devising a procedure allowing to cap- 
ture what do the "survived" theories have in common. In 
this way, even if we won't succeed in finding out a unique 
theory s e 5, we will be able nevertheless to find out a 
class of features, providing physical properties, that s has 
to possess so that it will not be falsified by the limited 
set of observations. 



A. Inverse problems and quantum dynamics 

In the last decades new ways of performing "observa- 
tions", computer simulations, have emerged in many ar- 
eas of scientific research and have provided new sources of 
inverse problems. In this work we focus on one in partic- 
ular, and namely on the estimation of spectral functions 
of many-body quantum systems starting from imaginary 
time correlation functions computed in Quantum Monte 
Carlo (QMC) simulations. Indeed, the study of dynam- 
ical properties requires the evaluation of spectral func- 
tions s(w): 
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A and B being given operators acting on the Hilbert 
space of the system whose Hamiltonian operator is H . 
The brackets indicate expectation value on the ground 
state or thermal average. Unfortunately, general meth- 
ods to obtain exact real time evolution are not known, 
thus forcing us to deduce a theoretical s{uj) from "obser- 
vations" of imaginary time correlation function /(r) = 
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{e^'^ Ae~^'^ 13) which can be computed by QMC meth- 
ods. The scenario is that of an inverse problem with the 
formal appearance of an integral equation 

/ + 00 
du:lC{T,uj)s{uj) (2) 
-oo 

where for example, at zero temperature, /C(r, oj) = 
6'(w)e~'^", 0{ijj) being the Heaviside distribution. The 
difficulty lies in the fact that the kernel /C(r, w) is a 
smoothing operator and QMC techniques are based on 
a discretization of the imaginary time domain, with time 
step i5t, allowing an estimation of /(r) only in corre- 
spondence with a finite number of imaginary time values 
{0,(Sr, 2(Sr, ...,Mr}, 

^={/o,/i,...,/z}. (3) 

In general F is obtained as an average of several QMC 
calculations of /(t), each affected by statistical noise and 
which are used to estimate the statistical uncertainties 
{cr^o , (T/j , . . . , CT/, } associated with J-. The task is then 
to evaluate s{u)) starting from limited and noisy data. 
Often sum rules provide useful help, either imposing ex- 
act constraints on s{u)) or allowing to perform additional 
QMC measurements which provide estimations for some 
moments of s(w): 

/ + 00 
dwa;"s(a;), n G Z}. (4) 
-oo 

For example co — {AB) may be easily estimated in equi- 
librium QMC simulations with an associated statistical 
uncertainty. Moreover some a priori knowledge may be 
assumed such as the support, non-negativity or some 
further properties. 

The task of facing the problem in equation ([2]) , namely 
an analytic continuation problem, has already been inves- 
tigated: the Maximum Entropy Method^ (MEM) is the 
most widely popular strategy developed; in the realm of 
bulk quantum fluids it has provided only qualitatively in- 
teresting results^'^. Another method has been proposed, 
the Average Spectrum Method^ (ASM), which has been 
recently applied to lattice spin models^ but also to realis- 
tic off-lattice systemaiS. The ASM strategy captures the 
idea of extracting physical properties using only observa- 
tions. Here we want to insist on this idea, introducing a 
new way to deal with the observations and their uncer- 
tainties strongly inspired by the falsification principle. 
As explained below, to do this we need a space of mod- 
els S, containing a wide collection of spectral functions 
consistent with any prior knowledge about s(a;), a falsi- 
fication procedure relying on the QMC "measurements" 
V = { J^, C} and a strategy to capture the accessible phys- 
ical properties of s{ijj). The strategy we are going to de- 
scribe in the following relies on Genetic Algorithms'^ to 
explore S and falsify its elements; for this reason we have 
called it the Genetic Inversion via Falsification of Theo- 
ries (GIFT) method. A preliminary application of this 



strategy to the determination of the dynamical structure 
factor S{q,uj) of liquid '^He was presented in RefjlS: here 
we explain in detail the GIFT method and we present 
several applications to the Helium system. 

The structure of the paper is the following. In section 

II we introduce the genetic inversion strategy; in section 

III we show applications of the GIFT method on several 
Helium systems; in section IV we present tests on the 
reliability of the GIFT method on known spectral models; 
section V contains our conclusions. 



II. GENETIC INVERSION STRATEGY 

A. The space of models 

In our mathematical framework S contains a wide 
class of step functions, providing a compromise between 
the possibility of suitably approximating any model of 
spectral function and the feasibility of numerical opera- 
tions inside it. In the typical case {A = B'') when s{lu) 
is known to be real-valued, non-negative and the zero- 
moment sum-rule holds, we rely on models s of the form: 

s{uj) differs from the physical spectral functions by a fac- 
tor Co, the zero-moment sum rule, which belongs to the 
set of observations and its role will become evident below. 
We introduce a discretization of the codomain: 

Sj & NU {0}, (6) 

to make the space finite, and we use the characteris- 
tic function XAj{^) of the intervals Aj = [LUj,LUj^i), 
{luq, ...,Ci;jv„} being a partition of width Aw of an inter- 
val of the real line larger than the hypothesized support 
of s{uj). M provides the maximum number of quanta of 
spectral weight available for the ensemble of the intervals 
A,. 

B. The falsification principle 

Once we have defined the space of model spectral func- 
tions, we have to devise a practical strategy to implement 
the falsification principle. We have to rely on the QMC 
estimations of the imaginary time correlation functions 
(and, in general, also of the moments) : 

J- ={/o, /!,...,/;}. (7) 

Such numbers are averages evaluated during a simula- 
tion and appear together with their estimated statistical 
uncertainties {cr/j, , crj^j , . . . , cr^j }. In typical approaches, 
such information are dealt with inside the framework 
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of Bayes' theorem; they provide the key ingredients to 
build up the a posteriori probabiUty^ to be maximized, 
together with some a priori probabihty, to extract the 
most probable spectral function. 

On the other hand, we find natural to suggest a 
new way of exploiting the information contained in 
{cT/o, CT/j , . . . , CT/, }: any set J^* equiv alent to i.e. such 
that I/* — fi\ is of the same order as cry., could be a re- 
sult of another hypothetical simulation. Falsifying the 
elements of S should require not only compatibility with 
J- but also compatibility with a vast population of J-* 
equivalent to the set ([7]) of data. At a practical level our 
aim is twofold: on one hand we need a recipe to generate 
equivalent sets J-* , on the other hand we have to use the 
generated F* to falsify the elements of 5. At the sim- 
plest level we have addressed the generation of the sets 
by sampling independent Gaussian distributions cen- 
tered on the original observations, with variances which 
correspond to the estimated statistical uncertainties so 
that a generic element T* is: 

^ {/o+e5, h+e\, ■ • ■ , fi+en = {/o*, ft. • • • , ft} (8) 

being random numbers sampled from Gaussian distri- 
butions with zero mean and variances equal to cr^ . . We 
stress that the very idea of exploiting the statistical un- 
certainties in the observations for generating equivalent 
sets J-* is the main difference with respect to preexisting 
statistical approaches to inverse problems. Note that in 
presence of more complete information in the observa- 
tions, like an estimation of the full covariance matrix, E, 
the generation of the equivalent sets F* can be readily 
generalized by sampling an {I -I- l)-variate normal distri- 
bution with the following probability density function: 
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exp [-\{F* - FyE-\F* - F)] 
(27r)^det(S)5 



(9) 



standard methods to perform efficiently this task are 
known (see for example Refflll). To be precise, being 
the set F a QMC estimation of the (unknowable) exact 
correlation function values, the probability density func- 
tion p{F*) in the previous equation does not represent 
the (unknowable) exact distribution of the observations; 
this situation is typical of the inverse problems we are 
considering and of the statistical methods developed in 
this context: we have to rely on what we have observed. 
Therefore, it is important to estimate the statistical un- 
certainties connected with the application of a particular 
statistical method. This can always be obtained com- 
paring models coming from the analysis of independent 
observations. 

The key point is then to falsify the elements of S rely- 
ing on each one of these sets: compatibility means small 
deviations from the observations. Thus, given the set F* , 
a very simple measure of the compatibility of a model 
with this set of observations can be obtained by comput- 
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The normalization of our models requires the multipli- 
cation of s(lo) by Co, the zero moment, which however 
belongs to the set of observations; a natural choice is 
to sample also its value analogously to how we treat F. 
This is the reason why a factor Cq appears in (fTOt . Each 
member F* of equivalent data leads to a different model; 
let us call Sfe the model found with the A;-th member F* . 
Each one of these models cannot be trusted to be the 
solution of the inverse problem, being at least partially 
biased by the particular F*; in other terms we can say 
that each one of these models will posses spurious in- 
formation, presumably different in each model, together 
with some physical information. An averaging procedure 
is therefore the simplest way to filter out the spurious in- 
formation and to reveal physical information, which con- 
sist in the common features among the models which have 
not been falsified. Therefore we take as the reconstructed 
spectral function the average 



SgiftH - 



k=l 



(11) 



where f\fr is the number of equivalent random set of F* 

(k) . ^ 

used in the computation and Cq is the Cq used in the 
k-th reconstruction. 

The average procedure in the definition of SciFTi^) 
points toward some similarities between our strategy and 
that of ASM. However, in the light of the falsification 
principle, the two methods are fairly different: in or- 
der to obtain its "solution" , ASM averages over spectral 
functions obtained by exploring model-space regions via 
a local Metropolis random walk based on a probability 
distribution;® thus ASM regards all the sampled spectral 
functions, also those with a low probability (of course, 
they are rarely sampled) as models which have not been 
falsified. Moreover, with ASM the statistical uncertain- 
ties in the observations play a different role, appearing 
only in the definition of the probability. Another issue is 
the algorithm used to explore the space of models; as ex- 
plained below, GIFT uses a non-local dynamics induced 
by a stochastic evolutionary process instead of a local 
Metropolis random walk which, in principle, could suf- 
fer from ergodicity problems, being the high probability 
model-space regions not guaranteed to be simply con- 
nected. Which of the two approaches, GIFT and ASM is 
superior might well depend on the specific inverse prob- 
lem or, in practice, it is possible that the reconstruc- 
tion abilities of GIFT and ASM are essentially equiva- 
lent. ASM has never been applied to the "^He case, it will 
be interesting in the future to compare the two methods 
on the same inverse problem. 

The natural scale of A(s) is provided by the value 
S — 7^ X]j=o ''^z ■ models s such that A(s) << i5 
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may provide unphysical overfitting. In our statistical 
approach to inverse problems there are two procedures 
which preserve from overfitting. The first one is that, 
given a set F* , the exploration of the space of models 
S should be stopped when a model s(aj) is found such 
that A(s) ~ (5; a further reduction of A will only repre- 
sents the intention to give to F* a strong belief, which 
is incompatible with the statistical treatment of the ob- 
servations in our strategy. The second procedure is even 
more relevant and in some sense it is intrinsic to our strat- 
egy: given an F* the reconstructed model s(w) contains 
some spurious information, but these information will be 
averaged out in Sgift{'-^)- 

At this point the following question urges an answer: 
How can we practically explore 5? We want to stress 
that the answer to this question is not the key point 
of our strategy; as we explain in the following, we have 
implemented genetic algorithms as efficient algorithms 
to explore our huge space of models, S. There could 
be inverse problems and different space of models which 
could be more efficiently explored with other algorithms; 
obviously, our strategy, which consists mainly in the new 
statistical treatment of observations, can be applied also 
in these cases. 



C. The fitness and the genetic dynamics 

Genetic algorithms (GA) provide an extremely efficient 
tool to explore a sample space by a non-local stochas- 
tic dynamics, via a survival-to-fitness evolutionary pro- 
cess mimicking the natural selection we observe in nat- 
ural world; such evolution aims toward "good" building 
blocks^'^ which, in our case, should recover information 
on physical spectral functions. The fitness of one partic- 
ular s{uj) should be based on the observations, i.e., on the 
noisy measured set V = {F,C}. But as explained before, 
taking into account the estimated statistical noise of V, 
any set V* compatible with V provides equivalent infor- 
mation to build a fitness function. Thus in our GA any 
random set V* — {J-*, C*}^"* gives rise to a fitness, which 
simply compares "predictions" of theories and "observa- 
tions" : 



*p*(s) = -A(s)-^7„ 



(12) 



In (|T2)) the free parameters 7„ > are adjusted in order 
to make the contributions to coming from J-* and 
from C* of the same order of magnitude: the idea is that 
we are not allowed to prefer some particular observation 
among the others, thus they should have the same weight 
in the fitness. If it happens that one c„ is exactly known, 
no error is added making c* = c„. We stress that 
provides the simplest and the most natural definition; 
moreover, as explained below, our GA uses <^t>* only to 
order models in ascending fitness, thus any alternative 
definition of which provides the same ordering will 
give rise to an identical GIFT algorithm. 



GA are well know procedures characterized by well de- 
fined (genetic-like) operations on populations of candi- 
date solution to optimization problems in applied mathe- 
matics. For basic nomenclature and standard implemen- 
tations one can refer to textbooks (e.g. see Il2i) . Here 
we simply sketch our particular realization related to the 
space of models we have defined. In our GA, we start 
randomly constructing a collection of s(w); each s(w) is 
coded by integers, Sj in equation ([5]). The genetic dy- 
namics then consists in a succession of generations during 
which an initial population, consisting of A/V individuals, 
is replaced with new ones in order to reach regions of S 
where high values of the fitness exist, for a given T)* . 

In practice, in the passage between two generations a 
succession of "biological-like" processes takes place, and 
namely selection, crossover and mutation. The selection 
procedure is meant to make individuals with large fit- 
ness preferentially be chosen to give rise to the next gen- 
eration: we achieved this by ordering the population in 
ascending fitness and selecting the fc-th individual with 
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where r is an uniform random number, r e [0,1), and 
[. ..] is the integer part; the non linear dependence of 
/c on r ensures that individuals with large fitness are 
preferentially selected. The crossover then operates on 
two selected s(w), the father and the mother, exchang- 
ing subparts of their total number of quanta of spectral 
weight, TW, to generate two sons. We have used a special 
single point crossover by sampling a random integer, w, 
between and and by exchanging w randomly cho- 
sen quanta of spectral weight between the father and the 
mother. In this way, the second equation in ([5]) is auto- 
matically satisfied, implying that the zero-moment sum 
rule is also satisfied. Each exchanged quantum remains 
in the original frequency bin as in its parents, thus en- 
suring that strong features present in both parents tend 
to persist in the sons. Successively, with a given proba- 
bility, mutation takes place on the two new individuals, 
i.e., a shift of a fraction of spectral weight between two 
intervals Aj. This is repeated till a new generation of 
s(cli) replaces the old one, with the exception of the s(w) 
with the highest fitness in the old generation which is 
cloned (elitism). The number of individuals in the new 
population is constantly reduced by about 5% at every 
generation till M-s is equal to a given minimal value; from 
this point over, the number of individuals N-s in the new 
generations is kept constant to this minimum value. The 
discarded individuals are those with the smallest fitness 
in the population. This is done to start the genetic evo- 
lution from a wide variety of possible models, but to not 
dissipate computational time on falsified spectral func- 
tions. 

In our context the GA dynamics performs the falsifi- 
cation procedure: only the s(w) with the highest fitness 
in the last generation provides a model, 5^(0;), which 
has not been falsified by V* . The evolutionary process 
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stops either after a reasonable maximum amount of gen- 
erations, A/g, or when a model s{u!) is found such that 
A(s(w)) ~ S. Many independent evolutionary processes 
are generated by sampling different V*, thus obtaining 
a set made of the elements Cq'^'' Sfc(a;); at this point, as 
explained above, an averaging procedure inside this set 
appears as the most natural way to extract physical in- 
formation. 



III. RESULTS FOR HELIUM SYSTEMS 

We are ready now to present applications of our ap- 
proach on physical systems. Long Monte Carlo runs have 
been performed in order to get imaginary time correla- 
tion functions with a typical statistical uncertainty of 
0.1-1%. For bulk superfluid *He most of the simulations 
have been for iV = 64 and N = 256 atoms moving in a 
cubic box, but also N = 128 and A'^ = 512 have been 
studied; for solid ^He the hep lattice with N = 180 and 
N ~ 448 lattice position have been used. Imaginary time 
correlation functions have been computed for instants 
Ti = ISt, I — 0,..,lmax = 60 in the superfluid phase 
and Imax = 30 in the solid phase, spaced by 5t = 1/160 
K~^. In the present applications on "^He systems the 
whole covariance matrix has not been computed, thus 
equivalent sets T* have been sampled simply by using 
the procedure in equation We have tested on a 

very simple quantum system that the knowledge of the 
whole covariance matrix for the imaginary time correla- 
tion function does not improve the reconstruction abili- 
ties of our strategyii^ when the kernel in equation ^ is 
of the form /C(t, w) = 9{uj)e~'^'^\ a similar conclusion was 
obtained in previous studies of superfluid ^He via MEM.^ 
All the results shown in this article have been obtained 
with the interatomic interaction of RefflTl. but some com- 
putations have been performed also with that in Reffisl 
as mentioned in the text. We have used the pair-product 
approximation^^ to express the imaginary time propaga- 
tor in the interval 5t — 1/160 K"-'^ which is known to be 
very accurate. For bulk superfluid '*He we chose 7n = 
Vn 7^ 1 (see equation (H^l), i.e., we have included only cq, 
which is the static structure factor, and the first moment 
sum-rule which is exactly known, ci — \q\^ /2m. For the 
extraction of the impurity branch and of the vacancy ex- 
citation spectrum we have only used the zero moment 
sum rule. Other parameters were fixed to Aw = 0.25 K, 
M = 5000, = 600 - 1600, initial value of 7V^ = 25000 
which is decreased down to the minimum value A/j — 400, 
as explained above; we have used about 10^ different sets 
V* and the number of generations for a given V* have 
been fixed to 10^. We have performed many tests with 
different choices of such parameters showing that none 
has a critical role once TV^^Aw is larger of the support of 
the reconstructed spectral functions. 
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Figure 1: (line) SGiFT{q,uj) for q = 0.783 A"^ and p = 
0.0218 A~^; (open circles) observed^^ dynamic structure fac- 
tor S{q, Lo) in liquid "^He for q = 0.7 A"^ at SVP and T = 1.3 
K. Notice the logarithmic scale. Notice also the difference 
between the wave vector of SGiFT{q,oj) and the one of the 
experimental available^^ dynamic structure factor; the exper- 
imental single particle peak position is known to increase by 
about 0.8 K in moving from q = 0.7 A~^ to g = 0.783 A"\ 



A. The dynamical structure factor of superfluid 
*He 

Our first case study is the determination of the dynam- 
ical structure factor S{q^ lo) of liquid bulk ''He. We have 
used the exact SPIGS method^'^-^^ to compute the inter- 
mediate scattering function F{q, r) at T = K near the 
equilibrium density, p — 0.0218 A^'^, and slightly above 
the freezing density, p — 0.0262 A""^; F{q,T) is simply 
/(t) when A = SMs chosen to be the Fourier transform 
of the local density operator A = p^ = X^il^i e""^'"' . We 
observe that our reconstructed Sgift{<1,^) exhibits an 
overall structure in good agreement with experimental 
data: a sharp quasi-particle peak and a shallow multi- 
phonon maximum are present (see Fig[T|) . Both features 
appear for the first time within an analytic continuation 
procedure applied to a QMC study of a many-body sys- 
tem in the continuum. Notice that it is not appropriate 
to compare the widths of the two sharp quasi-particle 
peaks in FiglT] in fact the experimental peak includes 
the broadening arising from instrumental resolution and 
the effect of the finite temperature; on the contrary, as 
explained in the following, the width of the reconstructed 
GIFT peak from a T = imaginary time correlation 
function is only a measure of the uncertainty in recon- 
structing its position. In FigI5]we show one SaiFTiq,^) 
in the roton region together with the excitation energies 
e{q) i.e., the position of the main peak as function of 
q. The uncertainties of e(g) correspond to the widths of 
the peaks a^: we have checked the consistency of such 
identification by performing independent QMC estima- 
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Figure 2: (a) and (b): Sgift{<1,'^) a-t q — 1.755 A ^ and 
p = 0.0218 A-^ (a) single quasi-particle (qp) peak; (b) mul- 
tiphonon (mp) contribution (notice change of scale). Lines 
corresponding to a SGiFTiq,^^) obtained with a nonzero en- 
tropic prior (»7 7^ 0) are also shown, (c) e{q) extracted at 
p — 0.0218 A~^ from the position of the qp (circles) peaks 
and the positions of the majcima of the mp contribution (tri- 
angles) are shown. The error-bars represent the 1/2-height 
widths, (d) e{q) and mp contribution extracted at p = 0.0262 
A~^. Lines in (c) and (d): experimental dat a'^^'^^ ; in the mp 
region in (c) the lower curve (dotted) represents the position 
of the maximum while the upper one (dashed) represents the 
1/2-height width. 



tions of F(q, r) and comparing the positions of the peaks 
obtained in SciFriq,^)', the distribution of the peaks 
displays a variance comparable to . 

In principle also a MEM-like algorithm could fit into 
the GIFT approach: it is enough to modify the fitness 
function by adding to $13* in equation (IT^ an entropic 
term —rjS^ with 



S — duj < s{uj) In 



s(w) 



m(cj) 



— s(uj) + m{u!) 



(14) 



S being the entropy as in Ref|6| and ij > a, free pa- 
rameter; m(cj) is the default model which in previous 
implementations^!^ has been chosen to be simply a con- 
stant in absence of any prior knowledge. This is not 
a faithful implementation of MEM because the entropic 
term is used in the context of the GIFT method and not 
within the framework of Bayes' theorem. Anyway, as we 
shall show, it provides results for the dynamical structure 
factor of superfiuid ^He very similar to the one appeared 
in literature.'^ The tests we are going to present in the 
following, in our opinion strongly suggest that the en- 
tropic term makes us loose a great deal of information, 
enforcing the smoothness of the spectral function as al- 
ready pointed out in RefS A good default model could 
improve these results, but such a priori knowledge gen- 
erally is not available. By using a constant as default 
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Figure 3: Evolution of the deviation l|10p during the stochas- 
tic evolution of the genetic algorithm for the reconstruction 
plotted in Fig[2ja,b) for 77 = averaged with respect to the 
sampled sets V* . The dashed horizontal line represents the 
value 5 = Ej=o '^/r 



model, m(uj), as in previous works^'^, for all wave vec- 
tors q we observed for the main peak of S{q, lo) a broad- 
ening (see Figl2]) strongly dependent on the choice of the 
parameter rj. This implementation of MEM provides re- 
sults comparable with those reported in references @ and 
0. The overall double-peak structure is typically lost in 
presence of the entropic term: no sharp peak is present 
but only a smooth S{q, uj) survives, as function of w, mak- 
ing the extracted excitation energy critically dependent 
on the value of ?/, thus introducing ambiguities in the in- 
terpretation of the results. In our original approach, i.e. 
without riS{s), we have checked that none of the param- 
eters (like A4, Aw, a, jn, •••) affects the class of features 
that we may trust to carry reliable physical information. 

As an example of the stochastic evolution, in Figl 
we show the deviation (jlOp as a function of the number 
of generations in the evolutionary process for the recon- 
struction plotted in Figl5Ja,b) for 77 = averaged on the 
sampled sets V*. One can see that the maximum num- 
ber of generations, Afg, we have used in this reconstruc- 
tion is optimal in reaching the "compatibility" condition, 
A(s) ~ 5 = Sj=o ' without overfitting. 

By integrating Sgift{<1,'^) we have access to quan- 
tities like the strength of the single quasi-particle peak, 
Z{q), and thus also to the contribution to the static struc- 
ture factor, S{q), coming from multiphonon excitations. 
Remarkably, Z(q) turns out to be in close agreement 
with experimental data (see upper Fig|4|), thus strongly 
suggesting that the shallow maximum in SciFTil,'^) at 
large energy carries indeed reliable physical information 
on the multiphonon branch of the spectrum. The po- 
sition of such multiphonon maximum (see FiglJJ:) is in 
qualitative agreement with experiments^^: as we show 
in the next section, within the present implementation of 
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Figure 4: (a) GIFT strength of the quasi particle peak Zici) 
as function of q at two densities and experimental data^^. (b) 
GIFT Static density response function at two densities 
and experimental data^^i^. Error bars of theoretical results 
are smaller than the symbol size. 

the GIFT method there is no possibility to recover the de- 
tailed shape of the spectral function like the multiphonon 
substructures given by high resolution measurements^^ of 
S{q,uj). In the lower panel of FigUwe show the static 
density response function xil) obtained evaluating the 
c_i from SaiFTil,^), the agreement with experiments 
is impressive, also near freezingSi.. 

The calculation of the excitation spectrum e{q) in su- 
perfluid '''He via QMC was addressed for instance in 
ReflU and in Refl27l. but here we are clearly much more 
ambitious because we aim to reconstruct the full spec- 
tral function. Our method is so powerful that it is able 
to reveal the effects of even fine details of the interatomic 
interaction. For example, the computed spectrum e{q) in 
the phonon region is about 0.7 K above the experimental 
value. We understand this as an effect of truncation of 
the inter-atomic interaction v{r) at a certain distance re- 
in most of our computations the interatomic potential is 
cut-off and displaced to zero at Tc = 6 A, and the equa- 
tion of state gives rise to an overestimation of the sound 
velocity by about 16%. We have performed some compu- 
tations with Tc — 14 A, in a simulation of iV = 512 ''He 
atoms and in this case the sound velocity turns out to 
be correct and now the theoretical s{q) at small q agrees 
with experiment within the resolution Alj. By comput- 
ing e{q) for many wave vectors in the roton region and 
averaging among the excitation energies nearby the ro- 
ton minimum, at p = 0.0218 A~^ we extract a roton 
energy of Er = 8.96 ± 0.47 K with the v{r) in Ref[l3 
and of Er = 8.67 ± 0.29 K with the v{r) in RefUl, a 
potential considered more accurate (Experimental roton 
energy^s at SVP Er = 8.608 ± 0.01 K). At p = 0.0262 
A~'^ we extract a roton energy of Er = 7.43 ± 0.34 K 
with the v{r) in Ref|l3 and of Er = 7.22 ± 0.27 K with 



Figure 5: Impurity ^He quasi-particle peak in superfluid 
''He at SVP for several wave vectors; in the inset the ex- 
tracted excitation energies are shown together with experi- 
mental data^S. 



the v{r) in Ref,18. (Experimental roton energy2£ at 24 
bar Er = 7.3 ± 0.02 K). 



B. Impurity and vacancy dynamics 

Another interesting test case is provided by liquid 
■^He in presence of one '^He impurity, in order to ex- 
tract the impurity branch which has been experimen- 
tally measured^S. Variational results for such branch are 
known^S but no results from exact QMC are available. 
This calculation requires the choice of A = e~^'''^'"^p, 
where fimp is the position of the impurity. In FigjS] 
we show the reconstructed spectral functions together 
with the estimated dispersion relation obtained from a 
simulation of V = 255 '^He atoms and one '^He atom 
at p = 0.0218 A^'^. The agreement with experimental 
data^^ is very good, thus providing a robust check of va- 
lidity of our approach. 

As a further application of GIFT we have studied the 
excitation spectrum of a single vacancy in hep solid ^He 
at p = 0.0293 A^"^, a density slightly above melting. The 
behavior of vacancies in solid ''He is of high interest be- 
cause vacancies and other defects are believed to have 
a key role in the possible supersolid phase of '^He at low 
temperature^ii^^. In order to apply GIFT to vacancy dy- 
namics the first step is the definition of a vector position 
Xy that allows to follow the "motion" of the vacancy in 
imaginary time during a SPIGS simulation. This prob- 
lem is much more difficult than the evaluation of the 
impurity branch, because the very definition of Xy is far 
from trivial due to the large zero-point motion of the 
atoms in the low density solid, x^ turns out to be a 
many-body variable, depending on all the vector posi- 
tions of ''He atoms, and even not free of ambiguities. We 
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Figure 6: Vacancy excitation spectrum in solid *He extracted 
from the SGiFT{Q,i^) of the vacancy-vector position Xv at 
p = 0.0293 A"^ in a hep lattice with N=447 particles along 
the principal symmetry directions: (a) FK, (b) FM and 
(c) FA. Two different algorithms have been used to obtain 
Xv: a coarse-grain algorithm"^"^ (CGR) and the Hungarian 
algorithm^i^ (HUN). Dotted lines represent the spectrum 
of a tight binding model (T-B) for the hep lattice22 obtained 
imposing the values of the band width along the FK and FA 
directions equal to the average between the values extracted 
from the two different algorithms. The arrow points out the 
vacancy-roton mode. 



have employed two different procedures to obtain : the 
coarse-graiE^2^ and the Hungarian^li^^. In Figl6] we show 
the vacancy excitation spectrum ey{q) extracted from 
the vacancy spectral functions {A ~ g-^?-^") obtained 
with the two methods. The results obtained with the 
two definitions of Xy are very similar, and at first sight 
make evident a picture of Bloch waves in the crystal; 
the agreement with a tight binding hopping model'^^ is 
good. Notice that £i,(g) represents the excitation energy 
with respect to the state with a vacancy with |g| = 0, 
i.e., ey{q) does not include the vacancy activation en- 
ergy. By fitting ey{q) with the tight binding expression 
we extract the vacancy effective mass in the different 
lattice directions: m^j^ = m^j^j = 0.46 ± 0.03m4 and 
TOp^ — 0.55 ± 0.1r7i4, where 1114 is the "'He mass; these 
values for m* are in agreement with the results obtained 
with a different method in ReflSTl 

The agreement of ey{q) with the tight binding model 
fails dramatically in the FM direction. In fact, at any re- 
ciprocal lattice vector the excitation energy should van- 
ish. On the contrary at the first reciprocal lattice vec- 
tor along FM the vacancy excitation spectrum does not 
vanish but it reveals a novel vacancy-roton mode with 
an energy of 2.6 ± 0.4 K and an effective mass of about 
= 0.46 m4. We have checked that this energy does 



not depend on the size of the system. Such behavior 
of £u(g) in the FM direction implies that the (non-zero) 
minimum is not a consequence of the lattice periodic- 
ity but it is related to correlated motion of particles like 
in superfluid ^He. It is interesting that neutron scatter- 
ing from hep ^He gives an unexpected excitation mode 
beyond the phonon modes exactly in the FM direction 
with a roton-like mode at the reciprocal wave vector'^'^. 
The experimental energy of such roton mode is about 
4.4 times larger than what we find; so it is unclear the 
connection between our mode and experimental data. A 
larger vacancy roton energy might arise in presence of 
clusters of vacancies. By analyzing the contributions to 
/(t) {e^'^Ae^^'^A'') with A = e"'"''^", one can see 
that the vacancy-roton mode is connected to motions 
of the vacancy between different basal planes. The fun- 
damental difference between in-basal-plane and inter- 
basal-plane motions is that the lattice position in the 
first case is a centre of inversion whereas this is not so in 
the second case. The fact that hep is not a Bravais lat- 
tice is fundamental in this respect. We have verified that 
in bcc crystal and in a two dimensional triangular lat- 
tice, both Bravais lattices, no such vacancy roton mode 
is present. 



IV. TESTS ON KNOWN SPECTRAL MODELS 

In this section we show several tests of application of 
the GIFT method on known analytical spectral models 
suitably discretized and "dirtied" with random noise to 
"simulate" actual data. It will appear evident what we 
have already pointed out in the introduction: only some 
features of the exact solution can be consistently repro- 
duced; we have no possibility to reconstruct exactly the 
shape of s(w); on the other hand, access is granted to 
the identification of the presence of peaks and to their 
positions, to some integral properties involving s{uj) and 
to its support. 

The most natural test for the reliability of the GIFT 
approach we have developed is provided by a systematic 
study of Laplace inversion problems whose analytical so- 
lution is known. Our idea is to focus our attention on 
model functions of the form: 

A/p — A/p 

linear combinations of Gaussians multiplied by 0{uj), the 
Heaviside distribution, resembling qualitatively the ex- 
perimental results for spectral functions in condensed 
matter physics at T = 0. We may perform several 
tests varying the parameters Afp, number of maxima, 
{/ii, /^TVp}, positions of the maxima, {ai, ...,a^f^}, 
widths of the peaks and {pi, ...,p_^^}, the areas under 
the peaks. The Laplace transform /(t) of may be 
expressed in terms of the standard complementary error 
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function: 



erfc(z) = —j= 



+00 



die 



(16) 



whose values are tabulated, in the following form 



i=i 



V2 



In order to simulate the output of a typical QMC cal- 
culation, we define the measured imaginary time data 
^ = {/o, ■••/;} as: 



f,=f{j5T)+s, 



(18) 



where /{jSr) is evaluated from (jl7p . and ej are random 
numbers, mimicking the error bars affecting QMC data, 
following Gaussian distributions with zero mean and vari- 
ances, a'^. , comparable with the ones typically occurring 
in our QMC results [a^Jfj in the range 0.1-4 %). fj 
play the role of the output of QMC simulation; GIFT 
falsification uses Mr random sets F* ~ {/q,.../*} de- 
fined by 



(19) 



£* being Gaussian random variables with zero mean and 
variances which here, to be coherent with the applications 
we have presented, we assume to be equal to a'^. . 

Our aim is to compare (IT5|) with the GIFT result we 
obtain pretending that our knowledge about the imagi- 
nary time correlation function is limited to the discretized 
and noisy data T in (jl8p . and to other available infor- 
mation about the moments c„ which, inside these tests 
on analytically solvable models, can be evaluated from 
(ITSI) ; we will neglect now the error bars affecting the val- 
ues of the c„. The parameters we have employed in our 
GIFT reconstructions are listed in Table ID Obviously, 
the choices of the interval of the frequency space, of the 
resolution Acj (which fixes N^^), of the discretization At 
and of the number of points in imaginary time / are cru- 
cial for a specific spectral function one is trying to recon- 
struct and should be chosen consistently with the con- 
sidered model; the other parameters are not crucial for a 
correct functioning of the GIFT method and they have 
been chosen in order to falsify a wide variety of mod- 
els leaving the computational cost of the algorithm at a 
reasonable level. 

Also in the reconstruction of known models of spec- 
tral functions one can compare GIFT results with those 
based on the strategy of MEM by adding in the fitness 
function an entropic term, as we did with the dynamical 
structure factors in superfluid *He from QMC imaginary 
time correlation functions. 



A. Single peak reconstruction 

The simplest test-case is provided by the attempt of re- 
constructing spectral functions displaying only one peak 




Figure 7: (Color online) Single peak reconstruction. Upper 
panel: Two noisy imaginary time correlation functions ob- 
tained via (|18p from /'"'(t) (open circles) and f^''\r) (x sym- 
bols), which are the Laplace transforms of Sa(io) (dotted line 
in the middle panel: fj, = 10 and a = 0.1) and st{uj) (dotted 
line in the lower panel: n — 10 and a = 1). The inset is a 
zoom on one imaginary time instant. Middle panel: Sa(tj) 
(dotted line) and reconstructed 5'g/ft(i^) (red line) using in 
the fitness $17* only the first moment ci (i.e. 7ri = Vn 7^ 1); 
green and blue lines represent MEM-like reconstructions with 
different values of the ri parameter in the fitness (see legend). 
Lower panel: Sb{uj) (dotted line) and reconstructed SciFTi^jj) 
using in the fitness <l>i>* only the first moment ci (i.e. 7n = 
Vn / 1). 



at a given point fi with a width a. The upper panel of 
Fig. [7]makes evident the difficulty of the inverse problem: 
two functions with the same parameter fia = = 10 but 
different values of the widths, respectively Ua = 0.1 and 
afc ~ 1.0, in imaginary time domain differ by about 0.5%, 
of the same order as the typical QMC error bars. It is 
clear then that the information about the width of the 
peak is always strongly obscured by the noise. However, 
from the middle and lower panel of Fig. [7] it is manifest 
that GIFT reconstruction, obtained using in the fitness 
only the first moment ci (i.e. 7„ = Vn 7^ 1), is 
able to capture with an high accuracy the position of the 
peak in both cases. We note that, despite the difficulty 
mentioned above, in this simple case of a single peak, 
even the widths are remarkably semi-quantitatively re- 
covered. This ability is evidently lost when the entropic 
term with a constant default model is switched on (see 
FigEl). 
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Figure 8: (Color online) Double peak reconstruction for well 
separated peaks. Upper panel: Noisy imaginary time correla- 
tion function obtained via H18|) from /(r) (dotted line) which 
is the Laplace transform of s(a)) (dotted line in the lower 
panel, see text for parameters). Lower panel: s{uj) (dotted 
line) and reconstructed Sgift{ijj) (red line) from /(r) using 
in the fitness <[>x>* only the first moment ci (i.e. "fn ~ 
Vn 7^ 1); green and blue lines represent MEM- like reconstruc- 
tions with different values of the rj parameter in the fitness 
(see legend). 





name of the parameter 


symbol 


value 


SM 


number of bins in frequency space 




600 


SM 


resolution in frequency space 


Auj 


0.25 


SM 


number of quanta of spectral weight 


M 


5 X 10^ 


CF 


discretization in imaginary time 


St 


1/160 


CF 


number of points in imaginary time 


I 


60 


GA 


number of generations 


Ng 


10* 


GA 


initial number of models 


Ns 


2.5 X 10* 


GA 


final number of models 


N- 


400 


GA number of new random sets generated 


Mr 


10^ 



Table I: Typical parameters used with the GIFT method re- 
lated to: the space of models (SM), the correlation function 
(CF) and the genetic algorithm (GA). 



B. Double peak reconstruction 



In order to get closer to realistic physical applications, 
we try to reconstruct also spectral functions displaying 
a double peak. Inside such a double peak reconstruc- 
tion, we may check also the estimation of the integrated 




Figure 9: Integral properties for double peak reconstruction. 
Upper panel: from the exact s(t<j) (dotted line) and 

/o(ti^) obtained with the reconstructed Sgift{'^) in Fig. |8] 
(case = 0). Lower panel: 7_i(aj) from the exact s{lo) (dotted 
line) and I-\(lo) obtained with the reconstructed Sgift{i^) 
in Fig. E 



spectral functions: 



loiu) = / diu'siuj') , 



du! 



(20) 



/o(a;) provides information about the spectral weight un- 
der the peaks in s(a;); in particular in the uj range be- 
tween the two peaks /o(w) gives the information from 
which we have derived the strength of the single quasi- 
particle peak, Z{q), in our GIFT study of superfluid ^He, 
as we will show in the following section. On the other 
hand, the asymptotic value of for large uj provides 

the key to estimate the static response function xil)- 

In Fig. |S]we show a reconstruction of a spectral func- 
tion s{uj) for two well separated peaks {pi — 0.5, p2 = 0.5, 
/ii = 10, /i2 = 21, ai = 0.1 and a2 = 2.0) using in the fit- 
ness only the first moment Ci (i.e. 7„ = Vn 7^ 1); 
this is the typical fitness used in our reconstruction of 
spectral functions of superfiuid ^He. The corresponding 
/o(w) and I-i{uj) are plotted in Fig. [9] compared with 
the analytic results from (fT5|) . We observe that no ap- 
preciable difference emerges, with respect to the exact 
results, as far as the determination of the positions of 
the peaks, of the areas under the peaks, and of the c_i 
moment (see Fig. [S]) are concerned: the accuracy is very 
good; on the other hand, the shape of the reconstructed 
s{uj) has not to be taken too seriously because it be- 
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Figure 10: Double peak reconstruction for overlapping peaks. 
Upper panel: Two noisy imaginary time correlation functions 
obtained via (fTH)l from /'"H-^) (open circles) and /''''(r) (x 
symbols), which are the Laplace transforms of Sa(i^) (dotted 
line in the middle panel) and Sb{oj) (dotted line in the lower 
panel); see text for parameters. The inset is a zoom on one 
imaginary time instant. Middle panel: Sa(cj) (dotted line) 
and reconstructed 5g/ft(i^) using in the fitness $-d* only 
the first moment ci (i.e. 7„ = Vn 7^ 1). Lower panel: Sb{ij) 
(dotted line) and reconstructed Sgift{(jj) using in the fitness 
^j}* only the first moment ci (i.e. 7„ = Vn 7^ 1). 



longs to the class of properties whose determination is 
obscured by statistical errors and discretization in imagi- 
nary time. MEM-like reconstructions are not even able to 
detect the presence of two peaks; moreover the position 
of the maximum of the reconstructed spectral function is 
dangerously 77-dependent (see Fig. |8]). 

In Fig. [To] we consider two different spectral functions 
Sa(cLi) and Sb(aj), characterized by two overlapping peaks, 
whose Laplace transforms, in imaginary time domain, are 
plotted in the upper panel {pia — 0.5, p2a = 0.5, fiia = 
10, ^J.2a = 15, aia = 0.1 and a2a = 4.0; pu = 0.5, 
P2b = 0.5, fj.ib = 10, ^26 = 15, aib = 1.0 and a2b = 4.0) 
As discussed previously, the small difference, comparable 
with the (pretended) error bars, rules out the possibility 
of a reconstruction of the actual shapes. Nevertheless, 
the GIFT method succeeds in finding out the positions 
of the peaks with good accuracy even in this case in which 
the overlap between the two peaks becomes significant. 



Figure 11: Multiple peak reconstruction. Upper panel: Noisy 
imaginary time correlation function obtained via (|18p from 
/(r) which is the Laplace transform of s{lj) (dotted line in 
the lower panel). Lower panel: s{uj) (dotted line) and recon- 
structed Sgift{i^) from /(r) using in the fitness ^-d* only 
the first moment ci (i.e. 7„ = Vn 7^ 1). 



C. Multiple peak reconstruction 

Finally, we devise the following test: we try to recon- 
struct a spectral function s{uj) {pi — 0.5, P2 = 0.1, 
P3 = 0.2, Pi = 0.2, ^ii = 10, fi2 = 21, ^3 = 27, 
/i4 = 35, ai = 0.1, a2 = 2, ~ 4, 04 = 6), display- 
ing a main peak and a broad contribution at higher w, 
made of a superposition of three Gaussians, resembling 
qualitatively the shape of the multiphononic contribution 
in the dynamical structure factor of superfluid ^He. We 
have tested our method using the usual fitness function 
<I>X)* with only the first moment ci included (i.e. 7^ = 
Wn^l): the results are plotted in Fig. [TTl In Fig. [Tithe 
integrated spectral functions are plotted; from the com- 
parison between the exact /o('^) and the one obtained 
from the reconstructed Sgift{^) one can observe that 
the spectral weights under the main peak and the broad 
contribution are well reproduced. Also the large uj limit 
of /_i(a;) is in good agreement with the exact value. 

One can also study the effect of the noise in /(r) in 
order to check the GIFT ability in recovering correct in- 
formation on the true s(w). In Fig. [T3| we show two 
Sgift{'^) reconstructed from a noisy /(r) with a^. 10 
times and 50 times greater than in the test shown in Fig. 
111! Only in the second case, which represents a situation 
of very high relative noise ((Te .//j in the range 5-200 
%), information on the correct spectral function is sensi- 
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Figure 12: Integral properties for multiple peak reconstruc- 
tion. Upper panel: /o(t<-') from the exact (dotted line) and 
/o(ai) obtained with the reconstructed Sgift{i^) in Fig. 1111 
Lower panel: /_i(ai) from the exact (dotted line) and I-i{uj) 
obtained with the reconstructed Sgift{ijj) in Fig. \TT\ 

bly lost. This test show the robustness of GIFT against 
noise in the observations, being able to recover correct 
information on s{lo) with a noise level up to one order of 
magnitude greater than what can be easily obtained in 
typical QMC calculations of imaginary-time correlation 
functions. 

It is possible also to use GIFT with a limited informa- 
tion on /(r), which corresponds as usual to /(r) values 
for a discrete set of imaginary times, but without any 
added noise. The result of such GIFT multi-peak re- 
construction is shown in the upper panel of Fig. 1141 By 
comparing this result with that shown in Fig. [Tl]it is pos- 
sible to see that the two Sgift{i-^) are very similar thus 
ruling out the necessity of more accurate observations of 
/(t) at discrete imaginary times in order to improve the 
GIFT performance. By maintaining the noise level in 
fj to zero, we have also tried to increase the amount of 
information by using I = 240 number of points in imag- 
inary time with 6t ~ 1/640; the result of such GIFT 
multi-peak reconstruction is shown in the middle panel 
of Fig. [TH No substantial improvement can be observed 
with respect to the previous case in spite of an increased 
computational cost of the GIFT algorithm. The compu- 
tational cost of GIFT is increased also by considering a 
wider space of model spectral functions. In our last test 
we tried a GIFT multi-peak reconstruction without noise 
with Acj = 0.1, the number of bins in frequency space 
Nuj ~ 1500 and the "quantization" of spectral weight 
Ai — 10"*. The result is shown in the lower panel of Fig. 
[Til here the noise in SciFTi^jj) is higher because due 



Figure 13: Panel (a): fj values obtained with a^^ — 0.01 
and used by GIFT to reconstruct s{lj) as shown in panel 
(b). Panel (c): fj values obtained with a^^ = 0.05 and used 
by GIFT to reconstruct s(a;) as shown in panel (d). Panel 
(b): exact s(aj) (dotted line) as in Fig. Illl and reconstructed 
<5'g/ft(i^) using the noisy observation of /(t) in panel (a). 
Panel (d): exact s(cj) (dotted line) as in Fig. [TT] and re- 
constructed SciFTiijj) using the noisy observation of /(r) in 
panel (c). 

to the computational cost of the GIFT algorithm with 
this parameters we have only averaged over Mr = 160 
random sets. Also in this case we found no substantial 
improvement in Sgift(j-^) as compared to the standard 
case. 

In the previous tests we have not explored different 
variants of GIFT as, for instance, a basis set different 
from step functions; one cannot exclude the possibility 
that by using different variants of this method more in- 
formation could be obtained. 



V. CONCLUSIONS 

We have built up a new strategy to attack inverse prob- 
lems which has been used to study dynamics in quantum 
many-body systems from QMC simulations; we have ob- 
tained very accurate results in the ^He case, in the liq- 
uid and in the solid phase, even in presence of disor- 
der, providing major improvements with respect to pre- 
vious studies appeared in literature on this case. We 
have stressed the important point that the problem we 
have faced belongs to the huge class of inverse problems, 
a deep topic also from a fundamental epistemological^ 
point of view. The basic idea of the falsification principle^ 
guided us in our particular implementation of analytic 
continuation, but, of course, may provide important im- 
plications in many fields of scientific research. Moreover, 
every problem emerging in Physics or applied Mathemat- 
ics that could be cast in the form of equation ([2]), inde- 
pendently of the specific kernel in equation can be 
attacked with GIFT. The method can be extended to in- 
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elude different kinds of constraints on the spectral func- 
tion or additional information like cross correlations be- 
tween the statistical noise of /(r) at different imaginary 
times. Many variants of GIFT can be devised depending 
on the problem, for instance a basis set different from 
step functions ([5]) can be used or non uniform discretiza- 
tion in presence of problems with multiple time scales, or 
distribution of noise that is not Gaussian. 



Figure 14: s{uj) (dotted line) as in Fig. [TT] Upper panel: 
reconstructed Sgift{i^) from the exact /(r) (i.e. without 
noise) assumed to be known for / — 60 number of points in 
imaginary time with St = 1/160, with Auj — 0.25 and using in 
the fitness "I>i5* only the first moment ci (i.e. 7„ = Vn / 1). 
Middle panel: reconstructed SoiFTiijj) from the exact /(r) 
(i.e. without noise) assumed to be known for I = 240 number 
of points in imaginary time with St — 1/640, with Acj = 0.25 
and using in the fitness "I'd* only the first moment ci (i.e. 
"fn ~ Vn 7^ 1). Lower panel: reconstructed Sgift{oj) from 
the exact /(t) (i.e. without noise) assumed to be known for 
I — 60 number of points in imaginary time with St — 1/160, 
with Aoj = 0.1 and using in the fitness only the first 

moment ci (i.e. 7n = Vn ^ 1). 
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